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Abstract 

Equal time spin-spin and pair field correlation functions are calculated for a 
two-chain Hubbard model using a density-matrix numerical renormalization 



^^ ■ group approach. At half-filling, the antiferromagnetic and pair field corre- 

j^ ■ lations both decay exponentially with the pair field having a much shorter 

g 

I ■ correlation length. This is consistent with a gapped spin-liquid ground state. 

^ ■ Below half-filling, the antiferromagnetic correlations become incommensurate 

O ■ 

and the spin gap persists. The pair field correlations appear to follow a power 

S^ . law decay which is similar to their non-interacting {7 = behavior. 

C3 ■ PACS Numbers: 74.20Hi, 75.10Lp, 71.10-hx 



Typeset using REVTgX 



Materials such as (VO)2P207 ^ and Sr2Cu406 0,^ contain weakly coupled arrays of 
metal-oxide-metal ladders. In the nominal insulating state, the metal ions have spin one-half 
and the oxygens mediate an antiferromagnetic superexchange coupling. Calculations for a 
two-chain antiferromagnetic model with an exchange interaction J along each chain and J' 
across each rung have shown that the insulating state can have a spin gap [^|-0]. A t-J, 
t'-J' model was introduced to describe the hopping of holes in the doped system. Lanczos 
calculations for this model with J' > t' found enhanced pairing correlations on a 2 x 8 ladder 
doped near one-quarter filling |Q. A mean-field analysis of this model using a Gutzwiller 
renormalization of the matrix elements found a spin gap in the insulating state which 
initially increases with doping. In addition, the doped system near half-filling was found to 
have a mean-field superconducting order parameter with a modified d-wave structure. 

Another representation of coupled chains is provided by the two-chain Hubbard model 
0. Here we report results obtained for the ground-state magnetic and pairing correlations 
of the two-chain Hubbard model using a density-matrix numerical renormalization group 
approach [|r^ to study lattices of up to 2 x 32 sites. We find that at half-filling, the 
dominant correlations are antiferromagnetic, but that these decay exponentially because 
of a spin gap, which we calculate directly. The singlet pair field correlations also decay 
exponentially but with a short correlation length of the order of the lattice spacing. The 
ground state at half-filling is thus a spin liquid. As the system is doped away from half-filling, 
the antiferromagnetic correlations show an incommensurate structure and a peak develops 
in the magnetic structure factor at a wave vector proportional to the filling. The spin gap 
decreases as the system is doped, but persists down to an occupation of 0.75 electrons per 
site or less. The singlet pair field correlations are enhanced and exhibit a power-law decay, 
but the form of the decay is approximately £~^, where d, is the separation distance. This is 
the same power-law dependence as for the noninteracting, t/ = 0, system. 

We picture the two-chain Hubbard model as a ladder standing along the y-axis with its 
rungs along the x-axis so that 



i,X(7 i^cr i,A 

Here c] ^^ creates an electron of spin a at rung j and side A = 1 (left) or 2 (right), the 
hopping along a chain is ty, the hopping between chains on a rung is t^, and U is the on-site 
Coulomb repulsion. 

The density matrix formulation of the renormalization group [|10| provides a controlled 
approximation for the calculation of ground-state energies and correlation functions. We 
have used the finite-system method with open boundaries and have studied 2 x 16, 2 x 24, 
and 2 X 32 lattices, keeping up to 400 states per block. Details of the calculations will be 
published elsewhere. The discarded density matrix weight (truncation error) varies from 
5.5 X 10^^ to less than 10^^ for the 2 x 32 results shown here. In addition, we have examined 
the convergence of measured quantities as a function of the number of states kept and 
ascertained that the symmetries of the Hamiltonian are preserved. The maximum errors in 
the quantities shown here are at most a few percent, and in most cases are much smaller. 
Since the method works best with chains that have open boundaries at the ends, all results 
here have open boundaries. 

We have calculated the equal-time spin and pair field correlation functions Sxxi{i,j) = 
(Ml^Mly), DUhj) = (A.iA^), and D,.(^,j) = (A,,A^) with 

At - rW^ - r^ r'^ (2) 

Here Suii,]) and 812(1, j) measure the spin-spin correlations along a chain and between 
the chains respectively, and Dxx{i,j) measures the singlet pair field correlations in which 
a singlet pair is added at rung j and removed at rung i. In addition, Dyx{i,j) measures 
the pair field correlations in which a singlet pair is added to rung j and removed from the 
right-hand chain between rungs i and i + 1. The relative phase of the pair wave function 
across the ith rung to along one chain from i to i + 1 is given by comparing the phase of 
Dxxihj) to Dyx{i,j). This turns out to be negative, corresponding to the mean field result 



obtained in Ref. |^. However, the non- interacting U = result at a filling (n) = 0.875 is 
also negative. 

Because of the open boundaries at the ends, the system is not translationally invariant. 
Therefore, the correlation functions are dependent on both i and j, whereas in the ther- 
modynamic limit or with periodic boundary conditions, they are a function of only \i — j\. 
We have found that the spin correlations at all fillings and the pairing correlations at half- 
filling are not strongly dependent on the placement of i and j, so we choose them to be as 
symmetric about the center of the lattice as possible, in order to minimize end effects. The 
pair correlations away from half-filling have strong variation with lattice placement, so in 
order to minimize these effects, we average over a number of i and j for a given |« — j|. We 
will then discuss the correlation functions as functions oi i = \i — j\ with these procedures 
implicit. As we shall see, there will be discernible boundary effects, but the lattices are long 
enough that one can still extract the general behavior. We also discuss the average filling, 
defined as {n) = l/NJ2i,xa ''^iM where A^ is the number of lattice sites. 

In order to understand the nature of the spin-spin correlations we have calculated the 
magnetic structure factor S{qx, Qy) by taking the fourier transform of Sxxi{i,j). Since there 
are two chains, Qx can take on only the values and vr. For the purposes of the fourier 
transform, we take Sx\'{i) = for £ larger than the lattice size. This introduces only a 
small error since the ^aa' {£) decays exponentially with i, and has decayed by a factor of at 
least 50 at the maximum lattice separation. The correlation function ^(vr, g^) is shown for 
four different fillings, (n) = 1, (n) = 0.9875, (n) = 0.875, and (n) = 0.75 in Fig. |l]. These 
fillings correspond to doping the half-filled system with 0, 2, 8, and 16 holes. At (n) = 1, 
there is a single strong peak at qy = tt. As the system is doped away from half-filling, the 
peak at qy = TT is strongly suppressed and ^(Tr, qy) peaks at qy = {n)-n. The residual peak at 
g^ = vr is present only for even numbers of hole pairs. Therefore, we believe that it is a finite 
size effect and will vanish in the thermodynamic limit. Since the correlations are strongly 
antiferromagnetic across the chain, S'(0, qy) (not shown) is small, flat and does not change 
much with the fllling. 



At half-filling, (n) = 1, the magnetic structure factor shown in Fig. |] has a Lorentzian 
line shape corresponding to an exponential decay of the spin-spin correlations. This is 
illustrated in Fig. |^, in which the logarithm of Sn{i) is plotted versus i for various values of 
U. The oscillatory deviations from linear behavior at large separation are due to the effects 
of the open end boundary conditions. From the slopes of the straight line segments, we have 
determined the correlation length ^ and this is plotted versus U/ty in the inset of Fig. ^ As 
U/ty increases, ^ decreases, saturating at a value of order 3 lattice spacings for t^ = ty. The 
interchain spin-spin correlations Si2{i) decay with the same correlation length as Sn{^)- 
At half-filling, the pair field correlations also decay exponentially, but with a correlation 
length of order the lattice spacing or smaller. Thus the pair field correlations are negligible 
at half-filling. 

We have also carried out calculations for U = 8ty at a cross-chain hopping t^ = 0.5ty and 
1.2ty. The correlation length decreases with increasing tx, ranging from 10.4 lattice spacings 
at tx = OMy, to 4.3 lattice spacings at t^ = ty and 2.9 lattice spacings at t^ = l.2ty. We 
expect the correlation length to diverge as t^ -^ because on a single chain there is no spin 
gap and the spin-spin correlation function decays as a power law [|n|]. The largest spin gap 
and thus the smallest correlation length should occur in the large t^/ty limit, in which the 
system can be described as decoupled spin singlets on the rungs of the ladder f^. 

In order to calculate the spin gap directly, one can calculate the energy difference between 
the lowest 3^ = state and the lowest S^ = I state. For the 2 x 32 system, the energy 
gap at half-filling is given by Agpin = £'(32, 32) — £'(33, 31) where E{N-^, N^) is the ground 
state energy of the system with A^^^ up electrons and N^ down electrons. The spin gap as a 
function of U at half-filling for tx = ty is shown in the main plot in Fig. ^ for a 2 x 32 lattice. 
We have examined the spin gap as a function of lattice size in order to confirm that the gap 
is present in the thermodynamic limit. For U = 8ty, Agpin = 0.151ty on a 2 x 16 lattice and 
Aspin = 0.132ty on a 2 X 32 lattice. At f/ = 0, the spin gap vanishes because there are a set 
of degenerate, half-occupied states on the Fermi surface. The gap peaks at approximately 
U = 8ty then decreases with increasing U. The gap decreases for large U because it varies 



as J = 4:ty/U in the large U limit. In this same limit, the correlation length saturates since 
it varies as the spin wave velocity divided by the gap. The inset plot shows the spin gap 
as a function of filling for U = 8ty. The gap decreases fairly rapidly with filling, but is still 
present at (n) = 0.75. 

When the Hubbard ladder is doped away from half-filling, ((n) = 0.875) one sees from 
Fig. |l] that the antiferromagnetic correlations develop an incommensurate peak. The pairing 
correlation functions -Dxx(^) and Dyx{i) are shown in Fig. || for t^ = ty and U = 8ty at 
fillings of (n) = 1 and (n) = 0.875. At (n) = 1, the pairing correlations decay quite strongly 
compared to those at (n) = 0.875. The relative sign of Dxx{(^) and Dyx{t) is negative at 
both fillings, consistent with a modified (i-wave structure. In order to determine the strength 
of the pairing correlations, one must consider their ^-dependence at large distances. For a 
one-dimensional system, we expect that any pairing correlation will at best decay as a power 
of £ and can in some cases decay exponentially, as we have seen for the half-filled system. 
For two chains, one can compare with the the non-interacting [/ = ladder, for which 

Dxx{l) = {l/2Ttiy [2 - cos{2kfi0)i) - cos{2kf{7T)i)] . (3) 

Here A;/(0) = cos~^{tx + Ai)/2 and kf^Tc) = cos^^{tx — /i)/2 are the Fermi wave vectors 
corresponding to the bonding and antibonding bands of the two coupled chains with /i the 
chemical potential. In order to examine the decay of the pair correlations, we have made 
the log-log plot of Dxx{i) shown in Fig. |^. We compare the U = 8ty results with the infinite 
system U = results given by Eq. (^. Our results for the interacting 2 x 16 and 2 x 24 
lattices are consistent with Fig. |^ and thus while end effects are present, it appears that the 
equal time pair field correlations of the interacting system decay approximately as £~^. 

In summary, we have found that the Hubbard model on two chains exhibits spin-liquid 
behavior at half-filling. Both spin-spin and pairing correlations decay exponentially with 
the pairing correlations having a much shorter correlation length. There is a spin gap which 
increases as the interaction U is turned on, peaks, and then becomes smaller. When the 
system is doped away from half-filling, the spin gap decreases with filling but persists down 
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to at least {n) = 0.75. The spin-spin correlations become incommensurate and the d-wave- 
like pair field correlation are enhanced. The pair field correlations appear to decay as a power 
law with exponent close to —2, similarly to those in the non-interacting U = system. 

We have recently received a preprint by Tsunetsugu et al. [|12| in which the spin gap 
and binding energy for two holes is calculated for the t-J model on two chains via Lanczos 
diagonalization. 

The authors thank N. Bulut, T.M. Rice, A. Sandvik, M. Vekic, E. Grannan, and R.T. 
Scalettar for useful discussions. R.M.N, and S.R.W. acknowledge support from the Office of 
Naval Research under grant No. N00014-91-J-1143 and D.J.S. acknowledges support from 
the National Science Foundation under grant DMR92-25027. The numerical calculations re- 
ported in this paper were performed at the San Diego Supercomputer Center. The computer 
code was written in C++. 
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FIGURES 
FIG. 1. The fourier transform S{'7T,qy) of the spin-spin correlation function Sxy{i). Here 

tx = ty and U = 8ty and the calculations were made on a 2 x 32 lattice. 

FIG. 2. A semilog plot of the spin-spin correlation function for (n) = 1 as a function of U/ty 
for tx = ty on a 2 X 32 lattice showing the exponential decay of these correlation functions. The 
inset shows the correlation length in units of the lattice spacing taken from the slopes of the lines. 

FIG. 3. The spin gap Agpin as a function of U/ty for t^ = ty calculated on a 2 x 32 lattice at 
(n) = 1.0. The inset shows Agpin plotted as a function of filling (n) for U/ty = 8. 

FIG. 4. The pair field correlation functions Dxx{^) = (^xjAj.) and Dyx{i) = (Aj^jA^) versus 
i = \i — j\ for (n) = 1.0 and (n) = 0.875. Here tx = ty and U = 8ty on a 2 x 32 lattice. 

FIG. 5. Log-log plot of Dxx{() = {^xAIj) versus / = \i-j\ at (n) = 0.875. The U = results 
are taken from Eq. (3) and the dashed line has slope —2. 
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